Metastable states in the planar 2d XY model and dissipation in superfluid flow 
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We use the Metropolis algorithm to study the stability of superfluid flow in a model system, namely 
the two-dimensional planar XY model. Flow properties are examined by studying the behaviour 
of the system in meta-stable "twisted" states. We demonstrate the stability of superfluidity in 
this model and we discuss the Meissner effect and velocity quantization. We also study the critical 
velocity and dissipation by vortex creation and rotational flow and their dependence on the geometry 
of the system. An expression for the average superfluid velocity as a function of time, Vs{t), is 
obtained and compared with experimental results. 
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I. INTRODUCTION 

Interest in superfluidity has never waned. Quite the 
contrary, since the experimental realization of trapped 
atomic Bose-Einstein condensatesi*^ (BEG) activity in 
this field has intensified dramatically and the recent ap- 
parent discovery of a supersolid phase in solid helium'^ 
will most likely increase the level of interest. In addition, 
there is intense interest in the low temperature properties 
of bosons where many questions remain concerning the 
phase diagrams of model systems and the phase transi- 
tions from the superfluid phase to various exotic phases. 

An important question concerning superfluids is that 
of dissipation and the critical velocity. Perhaps one of the 
first things one thinks of at the mention of superfluidity 
is "flow without friction" . Very early on. Landau^ real- 
ized that the superfluid can be treated as a dilute gas of 
noninteracting quasi-particle excitations (phonons) and 
obtained the famous and often cited Landau stability cri- 
terion for superfluids^: 
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The dispersion relation, w(fc), is the energy of a quasi- 
particle excitation of momentum hk in the fluid, and v^ 
is the Landau critical velocity above which dissipation 
sets in. Clearly, if the dispersion relation u}{k) is linear, at 
least for small fc, the Landau critical velocity is finite. On 
the other hand, v^ =Q for quadratic dependence oiuj{k) 
on k. When v^ ^ 0, dissipationless flow at finite velocity 
is possible and the superfluid is said to be "stable" . 

To measure the critical velocity of a superfluid such as 
liquid ^Helium at a temperature T < T\ = 2.18K, one 
might, for example, pull an object through the super- 
fluid: The velocity at which a drag force starts to act on 
this object can be thought of as a critical velocity. Such 
experiments have been done& and find a critical velocity 



which depends on the superfluid density, p^, but not on 
the geometry of the system. Furthermore, this critical 
velocity was found to be of the order of a few tens of me- 
ters/second which is in good agreement with the Landau 
prediction, Eq. (^, when applied to Helium. 

Another way to measure the critical velocity is to set 
the superfluid in motion in a tube (or through an ori- 
fice) and measure the superflow velocity at which dissi- 
pation sets in. This is perhaps more appealing physically 
than the first experiment because it involves a fiowing 
superfluid and because one may make connections with 
persistent currents observed in superconducting systems. 
Such flow experiments have also been performed^ and 
flnd that the critical velocity is of the order cm/s or less 
and does depend on the size of the orifice! In addition, it 
was found that the higher critical velocities are obtained 
with the smaller orifices which is perhaps counterintu- 
itive. These results are surprising when compared with 
the Landau prediction. It is thus clear experimentally 
that these "critical" velocities are not equivalent. 

The widespread application of the Landau criterion as 
a test of stability, Eq. (QJ , has been criticized^ because 
it applies to the case of an object moving in the fluid and 
not to the case of the fluid itself flowing through orifices 
or in a torus. Perhaps a reason for the ubiquity of the 
Landau criterion is that in many cases one can calculate, 
if only approximately, the dispersion relation uj{k). It is 
much more difficult to calculate, even numerically, the 
critical velocity observed in flow experiments. 

In this paper, we shall adopt the viewpoint of refer- 
ence^ and after a brief review, we shall present our sim- 
ulation. As usual, we have in mind here the two-fluid 
model of superfluids. In the superfluid phase, the total 
density p is the sum of two densities. 
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where p„ (ps) is the normal (super) fluid density. The 
superfluid component does not carry entropy and flows 



without dissipation. The total particle current is then 



j = PsVs + PnVn, 



(3) 



where w, (vn) is the velocity of the superfluid (normal) 
component. To obtain the expression for i/,, we write the 
superfluid wavefunction as 



V'(r) = 
where tpoi^ is real. Then 
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where p is the momentum operator. The superfluid ve- 
locity is proportional to the average gradient of the phase 
of the wavefunction: Clearly, if the phase does not change 
coherently throughout the volume of the system, z/, = 0. 
Following reference^, we consider two situations which 
demonstrate fundamental defining properties of the su- 
perfluid. 

(1) A torus containing liquid ^iJe at T > Tx is spun 
around its axis at very low angular velocit}ii&. Eventu- 
ally the liquid will come to equilibrium with the moving 
walls. Reducing T below T^, the liquid goes into its su- 
perfluid phase and the superfluid component is observed 
to come to rest and, to conserve angular momentum, the 
torus and normal fluid gain angular momentum. The ex- 
periment"'^^ demonstrates the analogue of the Meissner 
effect in superconductors. 

(2) Starting with the same setup as above, the torus 
is spun at high angular velocitjsifi. The temperature is 
then reduced below Tx and the torus brought to rest. 
Eventually the normal component will itself reach equi- 
librium with the walls and come to rest. It can then 
be verified that the angular momentum of the station- 
ary torus is non-zero: The superfluid component is still 
flowing and may continue to do so for a very long time. 
This is the phenomenon of persistent dissipationless flow 
similar to persistent currents in superconductors. Such 
experiments have been done by several groups, see for 
example referencesil^'-^. 

These two fundamental properties may be understood 
with the aid of the velocity quantization condition flrst 
proposed by Onsager— 
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where kq = h/m is the flux quantum, h is Planck's con- 
stant and m the particle mass. Clearly, the closed in- 
tegration path must enclose a "hole" in the system, ei- 
ther a vortex or a physical hole, otherwise the path can 
be shrunk continuously to a point and only ji ~ sur- 
vives. Another important property of superflow, which 
will come into play below, is that it is irrotational. 
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FIG. 1: Qualitative form of the free energy (arbitrary units) 
versus velocity in the presence of superfluidity for a 64 x 64 
system. When the velocity satisfies Eq. ||5J, the free energy 
has a local minimum (for the values of Vs shown) which ex- 
plains the metastability of persistent superfluid flow and the 
Meissner effect. As Vs increases, the corresponding local min- 
ima get shallower and eventually disappear rendering the flow 
unstable. 



Equation ^ demands that persistent flow take place 
at a well defined value ukq- This is not an equilibrium 
situation since the superfluid can reduce its free energy 
by coming to rest. The persistence of the flow, at least 
for Vs below a critical value, means that the fluid is in 
a metastable state which will eventually decay into the 
equilibrium stable state at restSti^. As we shall see, the 
lifetime of this metastable state depends on the velocity 
itself and on T. Figure ^ shows qualitatively the form 
of the free energy as a function of velocitjiSii^. This also 
makes clear that transitions between local minima should 
be possible and lead to the loss of superfluid momentum, 
in other words the begining of dissipation due to exces- 
sive velocity. Such phase slips were indeed observed ex- 
perimentally-i&. We now make precise what was meant 
by "low" and "high" velocities in the discussion of the 
two experiments above. The Meissner effect takes place 
when the initial velocity corresponds to less than half a 
quantum and the superfluid, seeking the nearest velocity 
satisfying Eq. ©, comes to rest, i.e. excludes all flux. 
When the initial velocity is high, i.e. larger than half 
a quantum, the superfluid will seek the nearest veloc- 
ity satisfying Eq. ^ and settle into the corresponding 
metastable state. 

It was suggestedi^iii that vortex formation is behind 
these transitions and although this mechanism is gener- 
ally accepted, the details are still not well understood 
and no numerical simulations have been done. 

In the following sections we shall address these ques- 
tions numerically for a model system. The paper is or- 
ganized as follows. In section II we briefly review the 



two dimensional XY model. In section III we present 
our results for the non-equilibrium simulations including 
velocity quantization, superfluid density, transitions and 
flux dissipation, scaling of lifetimes with the geometry of 
the system and comparisons with experiments. Conclu- 
sions are in section IV. 



where (ij) denotes nearest neighbour sites and f3 = 1/kT 
with k the Boltzmann constant. Since the superfluid 
component does not carry any entropy, the increase in 
free energy when the superfluid is in motion is due only 
to its kinetic energy. Then, with ps the superfluid parti- 
cle density and taking Vs purely in the x direction, one 
may write for small v^, 



II. MODEL 

Addressing these questions with Quantum Monte 
Carlo simulations of '^He or of model systems such as the 
bosonic Hubbard model, poses very difficult algorithmic 
problems. As discussed above, putting the superfluid in 
motion requires a gradient in the phase of the wavefunc- 
tion. This introduces a complex phase in the Boltzmann 
weight rendering the simulation extremely difficult. 

Instead, we shall use classical Monte Carlo to simu- 
late the two dimensional planar XY model which has 
been studied extensively in connection with two dimen- 
sional superfluidity. This model exhibits the well known 
Kosterlitz-Thouless {KT) transition at Tkt- For T > 
Tkt, there is a condensation of dissociated vortices, the 
phase is disordered and the correlation function decays 
exponentially. For T < Tkt it has a spin-wave phase 
characterized by tight binding of vortcx-antivortex pairs 
and power law correlation function. There is no sym- 
metry breaking in this phase and, therefore, no magne- 
tization. In the language of superfluids, the absence of 
magnetization is equivalent to the absence of EEC. In ad- 
dition, the low energy excitations, the spin waves, have 
a quadratic dispersion relation: uj{k) ~ k^. On the face 
of it, this might be taken to imply an unstable superfluid 
according to the Landau criterion. However, it must be 
recalled the these spin waves are thermodynamic in na- 
ture, not dynamic: This model has no intrinsic dynamics, 
unlike the three component XY model which has been 
shown to posses dynamic spin wave excitations with lin- 
ear dispersioiji^. 

One may justify using this model to study the dissi- 
pation properties of two dimensional superfluids by re- 
calling that experiments on "^He filmaiS, show that the 
transition to the superfluid phase is indeed in the KT 
universality class and that these films have finite criti- 
cal velocitjii2i2&. In addition, most explanations of dis- 
sipation in superfluid flow^iSSiSi^ are based on vortex 
formation which is certainly a defining property of the 
two-dimensional XY model. It is, therefore, reasonable 
to use this model to address the questions of surperfluid 
stability and dissipation. The quadratic dispersion of the 
low-lying excitations and the absence of BEC make these 
questions even more interesting. 

Consider, therefore, a two dimensional square lattice 
with N = Lx X Ly sites. The partition function of the 
XY model is then given by 
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F{vs) ~ Fq + -LxLyPs 



(9) 



from which immediately follow the expressions for the 
superfluid momentum density. 
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and the superfluid particle density 

1 1 d'^F{vs) 



Ps = 
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It is clear from Eq. © that the lowest velocity must 
correspond to n = 1 and consequently, for Vs to be small 
enough to justify Eq. ©, the system must be large. 
Equation Ullfl is often expressed in words by saying that 
the superfluid density is the curvature of the free energy 
as a function of the superfluid velocity. That this is not 
quite true is clear from Eq. 10) and Fig. ^ which together 
emphasize that ps is the curvature of the parabola pass- 
ing through the first few local minima of the free energy 
as a function of Vs- However, for brevity, we shall refer 
to this as the curvature of the free energy. 

Recalling that F = — lnZ//3, it is easy to apply 
Eqs. (|10lll|l to the XY model and obtain 
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where the notation {i,j):x means the sum is performed 
over nearest neighbours only in the x direction (since we 
took Vs to be in that direction) . Equation H12|) allows the 
determination of p^ in a nonequilibrium flow situation^ 
while Eq. H13|l gives ps also at equilibrium where fs = 0. 
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FIG. 2: Configuration witli Va ~ nno, n = 5 on a 50 x 50 
lattice. This configuration satisfies tfie velocity quantization 
condition and locally minimizes the free energy for n = 5 at 
T = 0. The direction of the arrows gives the value of Oi. 



FIG. 3: The flux quantum, n as a function of time in Monte 
Carlo Sweeps (MCS) for 100 x 100 system at T = 0.25 < 
Tkt- We see that n flows to the nearest integer value clearly 
exhibiting the Meissner effect and velocity quantization. 



III. NON-EQUILIBRIUM SIMULATIONS 
A. Metastability and superfluid density 

In this section we study some of the non-equilibrium 
properties of the two-dimensional XY model. In partic- 
ular, we are interested in the transitions among the local 
minima of Fig.^ in other words the onset of dissipation 
in the superflow and its relation to the critical velocity 
and thermodynamic stability of the superfluid. The tran- 
sitions between the local minima are driven by thermal 
fluctuations as the system attempts to minimize its free 
energy. Different simulation algorithms will lead to differ- 
ent lifetimes of the metastable states. However, the scal- 
ing of these lifetimes with the geometry of the system and 
the velocity of the superfluid will be the same. We chose 
to use the single spin flip local Metropolis algorithm. We 
shall see below that agreement with experiments is very 
good. 

Non-equilibrium simulations are performed by placing 
the system in initial conflgurations corresponding to flow 
at a chosen initial superfluid velocity. Equation jsj shows 
that the superfluid flows if there is a phase gradient, while 
Eq. places a condition on the allowed values of Vg. 
Figure 121 shows an initial conflguration corresponding to 
a flowing superfluid with flux tikq = 5ko- It is also clear 
that this conflguration corresponds to irrotational flow. 

We flrst show that superfluid flow in the two dimen- 
sional XY model is (meta) stable by demonstrating the 
existence of the Meissner effect and velocity quantization. 
As discussed in the introduction, if the initial n is less the 
1/2, a "stable" superffuid will satisfy velocity quantiza- 



tion by coming to rest thus exhibiting the Meissner effect. 
On the other hand, if n > 1/2, the system will evolve to 
the nearest integer value of n changing its velocity in the 
process. Figure |2I shows simulation results for four ini- 
tial values of the ffux quantum, n = 0.4, 0.6, 1.4, 1.6 for a 
100 X 100 system at T = 0.25 < Tkt- The behaviour of 
single conflgurations is shown for n = 0.4, 0.6 while for 
n = 1.4 and n = 1.6 we show averages over 1000 configu- 
rations (the dashed lines). The figure shows clearly that 
the Meissner effect and velocity quantization are both 
present in the two dimensional XY model below Tkt- 
This means that the free energy does indeed have the 
form depicted in Fig. ^ otherwise there would be no rea- 
son for the flux quantum to get stuck at integer values 
as shown in Fig. O 

That integer flux conflgurations are metastable, rather 
than stable, is shown in Fig. ^ where transitions are 
seen. It is also clear from this flgure that transitions 
do not always change the flux by the same amount. The 
change in the flux, and consequently dissipation and the 
flnal metastable or stable conflguration, depend on the 
starting value: The larger the initial flux, the larger 
the change and, consequently, the smaller the value of 
the flnal flux. The reason for this will be discussed be- 
low where we will study these transitions in more detail. 
Such transitions between local minima, also referred to 
as phase slips, have been observed experimentally^^. 

We now compare the superfluid density using equi- 
librium and nonequilibrium measurements. The equi- 
librium measurements were done, as usual, by using 
Eq. p3|) and random initial conflgurations with n — 0. 
The nonequilibrium measurements were performed by 
putting the system in an n = 1 initial conflguration and 




2.5x10 



5.0x10 
time (MCS) 



7.5x10 



1.0x10 



FIG. 4: The flux quantum n versus MCS for 100 x 100 system 
at r = 0.5. It is clear that the final value of n depends on 
the initial value: It is not necessary for the system always to 
tunnel to the n = state, nor does it have to tunnel from n 
to n — 1. The initial values for the flux are n — 6,8, 10, 12. 
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FIG. 5: Histogram of the superfluid momentum for a system 
of size 128 x 128. Upper panel: T = 0.9 and 45 simulations. 
Lower panel: T = 0.825 and 30 simulations. See discussion 
in text. 



using Eq. 1)1 0|l and l|12(l to measure the superfluid density. 
For r < 0.7 the n — I metastable state is so long lived 
that even after ten million sweeps the transition does not 
occur (for 128 x 128 system). For such long lived states 
the measurements are simply performed as if the system 
is in an equilibrium configuration. However, as Tj^t is 
approached, the lifetime of the metastable state becomes 



very short and more care must be taken. In this temper- 
ature range we measure the superfluid density as follows. 
The system is put in the n = 1 initial configuration and 
its evolution is followed until it makes the transition to 
the n — configuration, performing measurements ev- 
ery few sweeps. This is done many times for the same 
temperature (40 to 300 times); of course the transition 
time can be very different from one simulation to another 
at the same T (this will be discussed below) . The super- 
fluid momentum densitySSij PsVs, from all the runs for the 
same temperature is then histogrammed to decide if it is 
even meaningful to calculate an average of this quantity. 
The results for two temperatures are shown in Fig. [S] 
The relative heights of the peaks is not important since 
the height of the peak at pgVg = depends on how long 
the simulation is allowed to run after the transition has 
taken place. What is important is the height of the peak 
at nonzero momentum compared to its width and also 
compared to the height of the histogram in the transi- 
tion region between the two peaks. It is clear that for 
the upper panel of Fig. it is meaningless to calculate 
the superfluid momentum in the metastable state, while 
for the lower panel this quantity is well defined. 

With the help of such analysis we calculate ps in the 
nonequilibrium situation and compare with the equi- 
librium values in Fig. The agreement is excellent 
which demonstrates numerically the correctness of this 
approach. For T > 0.875, we can no longer extract ps 
this way: The initial superfiuid velocity, Vs = 2tt/Lx 
(in units of h/m and for L^ — 128) is greater than the 
critical velocity at these temperatures and the superfiuid 
kinetic energy is quickly dissipated. In order to get closer 
to Tkt, the initial Vs must be smaller which is impossi- 
ble for this system size since already n — 1. To achieve 
lower Vs, a larger system would be needed; this is an- 
other demonstration of the important interplay between 
system size and Vs- 

Note that the excellent agreement between the equilib- 
rium and non-equilibrium measurements of ps (which we 
have also verified with simulations at higher velocities) 
does not support Eq.(ll) of referenceS2i which claims a 
strong dependence of ps on Vs . 

An interesting aspect of the nonequilibrium measure- 
ments of Ps is the possibility of getting direct evidence 
for the universal jump conditionS^, 



Ps{Tkt) — —Tkt- 

TT 



(14) 



Inil a method based on higher order derivatives of the 
free energy was presented as a way for direct detection 
of this jump. In the present method, it is clear that 
when the system is close enough to Tkt the superfiuid 
momentum will simply vanish discontinuously. Already 
for L = 128, this happens for T just above 0.875 (the last 
nonequilibrium point shown in Fig. ^ , whereas Tkt ~ 
0.9. 
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FIG. 6: The superfluid density, ps, from equilibrium and 
non-equilibrium measurements versus T. The error bars are 
smaller than the symbols. 



B. Vortices and dissipation 

We now turn to the dissipation mechanism, i.e. the ex- 
citations which take the system from one local minimum 
to another. Onsagesi^ and Feynmanii argued that the 
creation of large vortices in the flowing superfluid was 
responsible for these dissipative transitions. Although 
the general features of this idea are widely accepted, dis- 
agreement remains on the details. Here we shall study 
and confirm numerically that, in the two dimensional 
XY model, transitions between the metastable persistent 
flow states proceeds via the formation of large vortex- 
antivortex pairs oriented orthogonally to the direction 
of flow. We will describe how it happens and study its 
dependence on the width of the system and the flux ve- 
locity. 

Our numerical simulations demonstrate that, for T < 
Tkt, thermal fluctuations take place in such a way that 
the large scale band structures. Fig. |21 are maintained 
due to the spin stiffness. As individual spins undergo 
thermal fluctuations, the bands of approximately paral- 
lel spins fluctuate as large scale elastic objects maintain- 
ing their large scale form and spanning the system (in 
the y direction by choice). Deformations of these elastic 
objects cost energy proportional to the curvature of the 
deformation and to the spin stiffness, i.e. ps. Eventually, 
two bands whose spins differ by 27r are deformed enough 
to touch. When this happens, the 2tt difference between 
these bands loses its meaning where they touch, but away 
from the contact region the difference still exists. This 
topological ambiguity is in fact the vortex-antivortex ex- 
citation which can trigger the transition. This situation 
is shown in Fig. \7\ which is a snapshot of the system 



depicted in Fig. El at T = 1/3 after 1100 Monte Carlo 
sweeps. The square in Fig. 13 shows where two bands of 
spin differing by 2tt have touched. A vortex is visible 
just below the square and an antivortex just above it; 
also visible is a flux line going down from the vortex and 
connecting it to the antivortex (we have periodic bound- 
ary conditions). In two dimensions, for T < Tkt, the 
energy of such excitations is of the order of Inr, where 
r » I/y is the separation of the vortex-antivortex pair. 

When such contact is made, it may be very difficult to 
break and dissipation of a flux quantum can then proceed 
by zipping together the two bands into a single one thus 
eliminating a quantum of flux. The reason for this is that 
for T < Tkt the vortices are confined: It is favorable for 
a vortex to find an antivortex and annihilate. In the sit- 
uation of Fig. Q the logarithmic confining potential will 
pull them together with the vortex moving down and the 
antivortex moving up to meet and annihilate zipping the 
the two bands in their wake. In the case of non-periodic 
boundary conditions in the y direction, the vortex and 
antivortex will each be attracted to its image thus mov- 
ing outward and crashing against the wall. The outcome, 
in both cases, is that the two bands are zipped together 
and a quantum of flux disappears, i.e. dissipates. 

Note that this mechanism, while involving vortices, is 
not quite the same as that discussed in reference^^. In2S, 
the vortex is assumed to nucleate at a singular point 
where the wavefunction vanishes, [ipl — > 0, which implies 
a vanishing density at that point. In the mechanism dis- 
cussed here, only the phase, modeled by the angles 9, 
fluctuates: The transitions here are caused by phase, not 
density, fluctuations. 

whereby a single band snaps into two parts each of 
which subsequently shrinks and melts away. While we 
have observed such events in our simulations, it appears 
that in the temperature and flux quantum ranges we 
studied that the zipper mechanism dominates. 

The role of the velocity is now straightforward to de- 
scribe qualitatively. The higher the velocity, the more 
closely packed the bands become and thus the higher 
the likelihood that more than one contact region be es- 
tablished. Consequently, several vortex-antivortex pairs 
may be created with the result that the transition can 
dissipate more than one flux quantum. This explains 
what is observed in Fig. 01 where it is seen that more flux 
quanta are dissipated if the initial velocity is higher. If 
the velocity is very high, so many vortices are created 
that the transition is essentially immediate. 

It is interesting to note that whereas the superflow is in 
general irrotational, this is no longer the case while the 
system is undergoing a transition. This is clearly seen 
in Fig. [3 where one observes sheer in the flow: The x 
flow velocity in the zones where the two bands have not 
yet touched is larger than that where contact has been 
established thus producing sheer at the interface. 

Before discussing the scaling of the lifetime of the 
metastable states, we flrst show that the escape time 
from a local minimum is well deflned and may be mea- 
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FIG. 7: The configuration in Fig. 13 after 1100 MCS at T = 
1/3. The square shows the region where two bands of spins 
diflering by 2-k have merged. A vortex (antivortex) can be 
seen below (above) the square. 



FIG. 8: The superfluid momentum versus time in Monte 
Carlo Steps for a 128 x 128 system at T = 0.83. The line is 
a fit of the form paVs = a(l — tanh(6(t — r<;sc))). The dashed 
line is a fit to the transition of another configuration which is 
not shown for clarity. 



sured accurately. Figure |S1 shows a typical evolution of 
the superfluid monientum as a function of time in Monte 
Carlo Steps (MCS). The solid line is a fit of the form 
PsVs — a(l — tanh(6(i — Tesc))) which allows us to mea- 
sure the lifetime, t^sc-, of the initial metastable state. The 
dashed line is a fit of the same form to the evolution of 
another configuration which we do not show. We shall 
return to this figure in the next subsection. 

As mentioned in the introduction, the critical velocity 
of superfluid helium passing through orifices decreases as 
the opening size is increased. Here we show that similar 
behavior is exhibited by the two dimensional XY model. 
In Fig.l^we show as a function of the width of the system, 
Ly, the average time (in MCS) to make a transition from 
the n to the n — 1 state for T — 0.25,0.5. We see that, 
for Ly > 32, T decreases as a power law as the system 
gets wider. This decrease in r with increasing Ly may be 
understood qualitatively with the help of the mechanism 
shown in Fig. [7| Increasing Ly makes it easier, at con- 
stant spin stiffness, to bend the spin bands by the needed 
amount to make them touch since the curvature is smaller 
the larger the Ly. However, a straightforward application 
of thermal activation arguments fails to give the observed 
behaviour. The energy of a vortex-antivortex pair in a 
the superfiuid flow field is given by^^ 



2 ^y 



E = 27rpslnLy - {2n) psU—^ 



(15) 



where the distance between the vortex and antivortex is 
taken to be Ly as is seen in Fig. |7| and the superfluid 
velocity is 2mT/Lx. This gives a transition ratei^ 



voLye 



-pE 



(16) 



where vq is the attempt rate. The average lifetime is then 



^ ^ ^27rp,/3-lg-/3p,n(2^)- 



Ly/L^ 



(17) 



This result which predicts an exponential decay of t with 
Ly does not agree with the observed numerical resulti 



.28 



For very narrow systems. 



Ly — 



8, 10, 16 in Fig. it is 



seen that r does not decrease with increasing Ly . This is 
because the system is becoming one dimensional in which 
case the spins are not stiff and therefore r — > 0. 

Figure ITUl shows the dependence of r on i^s = 2'KnK.Q. 
The metastable states {i.e. persistent superflow) are 
shorter lived for larger velocities, fs, which is in agree- 
ment with experiments. The exponential decay with in- 
creasing velocity appears in Ea. ll7l but with an exponent 
which is too large. For the larger values of n, the decay is 
no longer exponential because local minima become too 
shallow and also because the simultaneous dissipation of 
fiux quanta becomes more important. 



C. Comparison with experiments 

Decay of persistent currents has been studied both 
in bulk and films. The results of references^SiSSiSi were 
found to be consistent with time dependence of the su- 
perfluid velocity in the form 

Vs = A- Bint, (18) 

where A and B are empirical constants. In other words, 
the observed dissipation was very slow. However, refer- 
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FIG. 9: The average lifetime r (MCS) versus Ly for 8 < Ly < 
300 and L^ ~ 128. The transitions are from n — 7 to n = 6 
(up triangles), n = 6 to n = 5 (down triangles), n = 5 to n = 4 
(squares). For these cases T = 0.5. Circles: transition from 
n = 12 to n = 11 at T = 0.25. The averages are calculated 
over between 50 and 10^ realizations. The dashed line is a fit 
giving T = 5.7 x W^ /Ly'^ , the solid lines are to guide the eye. 



FIG. 10: The average lifetime of metastable states versus the 
flux quantum, n, for 128 x 128 (squares) and 128 x 130 (tri- 
angles) systems at T = 0.25. The solid line is r = 5.4 x 
10^^exp(-2n), the dashed line is r = 1.9 x 10^'^exp(-2.276n). 



encel^ observed both slow and fast dissipation, the for- 
mer well described by Eq. ((TH|l while the latter much bet- 
ter described by 



Ai 



(i + Bity 



(19) 



where Ai , Bi and r are empirical constants. Equation 
(fT5|l fails to describe the fast decays. 

In this subsection we shall compare the dissipation in 
the XY model with these experimental results. Figure|Hl 
shows, for one configuration, the superfluid momentum 
as a function of time clearly displaying the dissipation of 
kinetic energy. The behaviour in this figure follows nei- 
ther Ea. (|18|) nor H19|) and does not resemble figures 6 or 
7 of reference^ (see Figs. 1 131 and II 41 in this paper). How- 
ever, Fig.|H|shows the behaviour of only one configuration 
dissipating exactly one flux quantum whereas in the ex- 
periments one presumably observes an average of such 
processes. In other words, the experimental situation 
represents many different regions of quantized flux mak- 
ing transitions at different times. What is observed then 
is the average dissipation as a function of time. To test 
this idea, we performed simulations of the type shown 
in Fig.|Slbut averaging over many configurations. Three 
such averages are shown by the points in Fig. 1111 This fig- 
ure strongly resembles Fig. 6 of referencei^ and the data 
are reasonably well described by Eq. (|19|l (not shown in 
the figure, the curves shown will be discussed below). We 
therefore see that, at least for fast dissipation, the XY 
model is in very good qualitative agreement with exper- 
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FIG. 11: Superfluid velocity versus time averaged over several 
realizations for a 128 x 128 system. Lowest curve: T = 0.92, 
100 realizations, middle curve: T = 0.875, 30 realizations and 
upper curve: T = 0.85, 30 realizations. See text after Eq. 12211 
for a discussion of the lines. 



iments adding support to the phase zipping dissipation 
mechanism discussed above. Slow dissipation is harder to 
study numerically because for long escape time it takes 
many more realizations to get good statistics. However, 
the simulations we performed for this case are also con- 
sistent with the experiments. 

To model the numerically observed time dependence of 



Vs displayed in Fig. ^2 one must examine the distribution 
of escape times, Tf,scj at a given T and Vg. One such 
distribution is shown in Fig. El The soUd curve is a fit 
of the form 



P(Te 



^T^^, 



-7''"e. 



where A, a and 7 are fitting parameters, 
of Fig. [El a = 3.7 and 7 = 3.1 x lO""*. 
superfluid velocity, Vs, is then given by 



Vsit) 



dTescVs{t)P{Te. 



(20) 

For the case 
The average 



(21) 



where P{Tesc) is normalized and where Vs(t) may be 
taken of the form discussed in Fig. |H1 However, to 
simplify the discussion and allow exact integration of 
Eq. (|^ . we may take Vs{t) = Vs{0)Q{Tesc — t) where 
Q is the Heaviside function. Equation 1211) then reduces 
to 



Vsit) 



A I dTe.,<,,e-^-- 



yl'r(l + a,7i) 



(22) 



where r(l + a,jt) is the incomplete Gamma function, 
A', a and 7 are fitting parameters. The solid curves in 
Fig. ^2 are fits to the numerical data using Eq. H22|) with 
(a = 0.08,7 = 2.77 x 10^^) for the lowest curve, (a = 
0.04,7 = 10~^) for the middle curve and {a = 0.02,7 = 
3.5 X 10^^) for the top curve. We see that agreement 
with the numerical results is excellent. As a further test 
of these ideas, we show in Fig.^la fit to the experimental 
data in Fig. 6 of referencei^ using Eq. l(^ . The curve 
agrees remarkably well with the experimental data, in 
fact much better than Eq. (|19|) , see Fig. 6 in referencai^. 

The escape time distribution. Eg. 1201 used in this anal- 
ysis is reasonable for fast to medium dissipation as can 
be seen in Fig.^] However, it is not clear that the same 
distribution gives a reasonable description for extremely 
long lived persistent flows of the type well modeled by 
Eq. (|18|) . In other words, the question is whether for 
very slow decays Eq. H22|) behaves like the experimental 
results, i.e. roughly linearly in Int. To verify this, we 
show in Fig. ^^ fits of Eq. (|^ to two data sets from 
figure 7 of referencei^. We see that even for the slowest 
dissipation reported, Eq. (j^ gives very good agreement 
with experiments. For the slow decay, though, it was 
necessary to take a < to get a good fit. This might 
cause concern in view of Eq. H20|) which would then di- 
verge for Tesc = 0. However, negative a is needed only 
for the very slow decays where t^sc is never zero. In fact 
the distribution of escape times for slow dissipation will 
most likely have a lower cutoff below which no dissipation 
is observed. Note from Ea. H20|) that 7^^ is a measure of 
the lifetime of the metastable states. The values obtained 
from the fits in Fig. ^J are 7""'^ w 2sec for the fast decay 
and 7~^ Ri lO^sec for the slow one. 

This discussion supports the view that the experimen- 
tally observed dissipation curves are averages over many 
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FIG. 12: The distribution of escape times from n — 6 to n — 5 
for a 128 x 128 system at T = 0.5. The sohd curve is a fit of 
the form Ea.l^ 
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FIG. 13: Superfluid velocity versus time. The circles are the 
experimental data in figure 6 of referencei-. The solid curve 
is a fit using Eq. H22|l . The fitting parameters are (q = 0, 7 = 
0.043). 



events like the one shown in Fig. IHlwith an escape time 
distribution similar to Fig. 1121 



IV. CONCLUSIONS 

We have presented a first study of the two dimensional 
planar XY model in topologically non-trivial metastable 
states. These twisted states represent a two dimen- 
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FIG. 14: Superfluid velocity versus time for films with differ- 
ent thickness at T = 1.45K. The symbols are experimental 
data from figure 7 of reference—. The solid curves are fits 
using Eq. ^. For 6.4 layers (a = 0.3,7 = 0.5) and for 9 
layers {a = -0.73, 7 = 0.001) 



sional superfluid under non-equilibrium persistent flow 
conditions which satisfy the superfluid velocity quantiza- 
tion condition, Eq. ((HJ. Transitions among these states 
change the topological quantum number, the quantized 
flux, and represent the dissipation of flux quanta in a 
flowing superfluid which is related to its critical velocity. 
The properties of these transitions were studied in par- 
ticular their dependence on the geometry of the system 
and the superfluid velocity. We showed that, in agree- 
ment with experimental results, the metastable states 
have shorter lifetimes and, therefore, lower critical ve- 
locities when the width of the system increases. Also in 
agreement with experiments, we showed that dissipation 
is more rapid the higher the initial superfluid velocity. 



The dissipation mechanism was identified and studied: 
When two bands of spin differing by 2tt deform and touch 
due to thermal fluctuations, a bound vortex-antivortex 
pair of length Ly is created. The confining force pulls 
the vortex and antivortex towards each other thus zip- 
ping together the two touching bands of spin into a single 
band. 

Using this mechanism and a functional form for the 
escape time distribution motivated by the numerical re- 
sults, we calculated the average superfluid velocity, Vs{t), 
and showed it to be in excellent agreement both with our 
numerical simulations. Fig. 1111 and with the experimen- 
tal results of referencei^, Figs. Upland [TH This provides 
support for the view that the dissipation observed experi- 
mentally is an average over several regions with, possibly, 
different velocities where dissipation of flux quanta takes 
place independently and at different times. 

It is interesting that the two dimensional planar XY 
model, with no condensate and no linear dispersion for 
the low-lying excitations still exhibits (meta)stable su- 
perfluid flow and agrees so well with experimental results. 

Finally we mention that the free energy landscape 
may be studied directly using, for example, the Wang- 
Landau^ algorithm. Integer windings in the two dimen- 
sional XY model were addressed in'^^. Very recentlji^, 
different topological sectors and the helicity modulus 
were studied for the four dimensional compact U{1) lat- 
tice gauge theory. 
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